Non-equilibrium relaxation and critical aging for driven Ising lattice gases 
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We employ Monte Carlo simulations to study the non-equilibrium relaxation of driven Ising lattice 
gases in two dimensions. Whereas the temporal scaling of the density auto-correlation function in 
the non-equilibrium steady state does not allow a precise measurement of the critical exponents, 
these can be accurately determined from the aging scaling of the two-time auto-correlations and the 
order parameter evolution following a quench to the critical point. We obtain excellent agreement 
with renormalization group predictions based on the standard Langevin representation of driven 
Ising lattice gases. 
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Driven diffusive systems represent paradigmatic mod- 
els that display non-trivial non-equilibrium stationary 
states, and hence constitute crucial test cases for vari- 
ous analytical and numerical approaches [EI3|. Already 
driven lattice gases with mere site exclusion and periodic 
boundary conditions yield anisotropically scale-invariant 
steady states. In one dimension, the asymptotic scaling 
properties of these asymmetric exclusion processes are 
governed by fluctuation-controlled exponents that differ 
from mean-field predictions; they belong to the noisy 
Burgers or one-dimensional Kardar-Parisi-Zhang (KPZ) 
universality class 0, Hj] . Related models have found var- 
ied applications in the theoretical description of biologi- 
cal systems [H ■ Adding nearest- neighbor attractive Ising 
interactions, as proposed by Katz, Lebowitz, and Spohn 
(KLS) [H, 0] to describe ionic conductors, promotes clus- 
tering of empty and filled sites, inducing a continuous 
non-equilibrium phase transition between a disordered 
and a striped phase in d > 2 dimensions. In the absence 
of the drive, one recovers the equilibrium Ising model 
ferromagnetic phase transition. Yet the correct long- 
wavelength description of the KLS Ising driven lattice gas 
at the critical point has been subject to a lasting debate 
Ambiguities in determining the critical exponents 
are in part caused by exceedingly slow crossovers to the 
asymptotic regime [13l |: also, it is crucial to implement 
the proper anisotropic lattice scaling, and infer scaling 
exponents sufficiently close to the critical point flij ]. 

Intriguing universal scaling features moreover arise in 
the non-equilibrium relaxation regime from an initial con- 
figuration that drastically differs from the final stationary 
state [lj|. When the order parameter field is conserved, 
one may relate the ensuing aging scaling exponents to 
those describing the non-equilibrium steady state flij ]. 
Investigating the aging and initial-slip scaling regimes 
thus provides an independent and often more easily ac- 
cessible means to determine the asymptotic critical expo- 
nents [17j] . We recently utilized the aging relaxation scal- 



ensuing aging scaling were also analyzed for the continu- 
ous phase transition in the contact process [I(J[2lj]. Here, 
we employ Monte Carlo simulations to study the non- 
equilibrium relaxation properties of the KLS model fol- 
lowing a sudden quench from the high-temperature dis- 
ordered phase to the critical point. We utilize the two- 
time density auto-correlation function and the temporal 
evolution of the order parameter in the aging regime of 
large systems to cleanly measure the associated critical 
exponents in two dimensions, see also Refs. 22, 23|. 

The KLS model consists of N hard-core particles on 
a Lii x L^ 1 lattice with periodic boundary conditions; 
"||" and "_L" denote directions parallel and perpendicular 
to the drive. The occupation number ni at any lattice 
site i is restricted to or 1; multiple site occupations 
are prohibited. In addition to this on-site exclusion, we 
consider an attractive nearest-neighbor Ising interaction 
H lnt — — ^S(ij) n i n h J > 0- The transition rate evolv- 
ing the system from configuration X to Y then is 



R(X ->• Y) oc e 



-/3[H iat (Y)-H iat (X)+eE] 



(1) 



ing to accurately confirm the exponent values for asym- 
metric exclusion process es 11811 ; for a related study of the 
KPZ equation, see Ref. [l9(. Critical quenches and the 



with inverse temperature /?; I = { — 1,0, 1} indicates par- 
ticle hops against, transverse to, and along with the drive. 
The bias E and periodic boundary conditions generate a 
non-zero particle current through the system; since the 
rates ((TJ) do not satisfy detailed balance, the system will 
eventually settle in a non-equilibrium steady state. For 
E = 0, and with rii = h(ai + 1), the KLS lattice gas re- 
duces to the equilibrium Ising model with spin variables 
(7i = ±1 and conserved Kawasaki kinetics. We employ 
the standard Metropolis algorithm with rates ((T|), setting 
J = 1 and limiting our study to the case E — ¥ oo, imply- 
ing that particle jumps with I = — 1 are strictly forbid- 
den. Time is measured in Monte Carlo sweeps, wherein 
on average one update per particle is attempted. 

The theoretical treatment of the long-time and large- 
scale properties of the KLS model uses a continuum de- 
scription in terms of a density field p(x, t) or equivalently 
the conserved spin density 4>(x, t) = 2p(x, t) — 1. Janssen 
and Schmittmann J24| as well as Leung and Cardy [25l | 
constructed a mesoscopic Langevin equation for the crit- 
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ical KLS model near the critical point, 



+ A 



(2) 
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and analyzed it by means of perturbativc ficld-thcorctic 
dynamic renormalization group (RG) methods near the 
upper critical dimension d c = 5. In the ordered phase, 
the system forms stripes oriented parallel to the drive, 
and correspondingly the critical control parameter is 
t = (T — Tq)/Tc —> 0, whereas c > 0. £ represents a 
coarse-grained driving field, with £ — > const as E — > oo. 
Since only the transverse sector becomes critical, just the 
associated Gaussian noise needs to be accounted for, with 
(£,±(x,t)) = and correlator (£j_ i(a?, t) £j_ t j(x', t')) = 
n±5(x — x')5(t — t')Sij = 1, . . . , d — 1). 

We are primarily interested in the density-density cor- 
relation function C(x, x' , t, s) = {p(x, t)p(x' , s))—p 2 , with 
the mean density p = (p(x 7 t)) = \ in our simulations. 
Assuming spatial and temporal translational invariance 
in the non-equilibrium steady state, one arrives at the 
general scaling form (with positive scale parameter b) [l| 

C(xj_,xpt,T,Lx,Ln) = (3) 
= b d+A ~ 2+ v C(x ± b, Xll b 1+A ,tb z , rb- l ' v , L ± b, ij|fe 1+A 



where t, xy, and x± now represent time and position dif- 
ferences. Here, v and z denote the transverse correlation 
length [2l| and dynamic critical exponents, the latter as- 
sociated with critical slowing-down, while 77 characterizes 
the power law correlations at criticality. The external 
drive induces spatially anisotropic scaling, captured by 
a nonzero anisotropy exponent A (in the mean-field ap- 
proximation A = 1), which in turn results in distinct 
exponents in the direction of the drive, v\\ = v{l + A), 
whence z\\ = z/(l + A), and r]\\ = (r) + 2A)/(1 + A) 
[H . In the RG analysis near the upper critical dimension 
d c = 5 as determined by the coupling associated with the 
drive £ , the static nonlincarity u becomes (dangerously) 
irrelevant, and the overall conservation law together with 
Galilean invariance fix the scaling exponents to all orders 
in the e = 5 — d expansion [3, Ei| ■ The resulting exact 
numerical values (henceforth referred to as " JSLC" ) are 
listed in Table U for e > 0, and specifically for d = 2. 

Early Monte Carlo simulation data [13-01 °f the criti- 
cal KLS model with infinite drive indicated discrepancies 
with these JSLC predictions. Leung subsequently ap- 
plied careful anisotropic finite-size scaling to his Monte 
Carlo data, and obtained agreement with the JSLC the- 
ory, also supported by later work 12|, |l3|. A debate 



on the validity and implications of his results followed 
14, [loj], and it was proposed 31-33J that the critical 
KLS model at infinite drive might not be adequately de- 
scribed by Eq. in the asymptotic limit, but instead 



TABLE I: Critical exponents for the KLS and RDLG models 
(precise definitions are provided in the text). 





JSLC - exact 


RDLG - 0(e 2 ) 




e = 5 - d 


d = 2 


e = 3 — d 


d = 2 




l + e/3 


2 


1 - 2e 2 /243 


0.992 


z 


4 


4 


4 - 4e 2 /243 


3.984 


V 


1/2 


1/2 




0.626 


n 








4e 2 /243 


0.016 


P 


1/2 


1/2 


Kl + A-e + r,) 


0.315 


C 


1 - e/6 


1/2 


(2- £ + ifl)/(4-!fl) 


0.253 



by the Langevin equation for the randomly driven Ising 
lattice gas (RDLG, or two-temperature model B) (34^ . 



d t <t>(x,t) 



+ A 



-V||£||(£,i)-Vx-a(z,*) (4) 
Vl) (j}(x,t) +uV 2 L (j)(x,t)^ 



rV 2 



Here the drive term ~ £ that is the origin of the non- 
vanishing particle current has been dropped, with the 
reasoning that the presence of anisotropic noise with dif- 
ferent strengths n± and n\\ is supposedly more relevant in 
the infinite drive limit. The fluctuations for the Langevin 
equation ([4]) are controlled by the static nonlinearity u 
with upper critical dimension d' c = 3. Dynamic RG 
methods had earlier been employed to determine the as- 
sociated critical exponents to two-loop order in a dimen- 
sional e = 3 — d expansion [3~4T - [36j ; the resulting explicit 
values are also listed in Table HI However, numerical data 
for various finite-size scaling functions by Caracciolo et 
al. favor the JSLC rather than RDLG description [33 |. 

We have performed extensive Metropolis Monte Carlo 
simulations for the critical KLS model, based on the rates 
(TTJ) to numerically determine the steady-state scaling ex- 
ponents. To account for the anisotropic scaling laws, 
the simulations should be run with properly scaled sys- 
tem extensions, where Lu = AL^ A , where A = A co- 
incides with the correct anisotropy exponent, such that 
A = const. If A 7^ A is chosen, A does not remain con- 
stant and will enter the finite-size scaling functions as 
an additional relevant variable [1, H3 > see a ^ so -R- ei - @] ■ 
Setting b = L^ 1 ^ 1+A \ t = 0, and letting Xj_,x\\ — > 0, 
Eq. © reduces to the steady-state finite-size scaling form 
for the density auto-correlation function at criticality, 



C(i,L||) 
Following Rcf. 



j- -(d+A-2+7j)/(l+A) n 



FS 



t/Ll 



(5) 



T, 



we determined the finite-size "critical 
± for three different system exten- 



temperatures 

sions (L\\,L±) = (54,24), (128,32), and (250,40) (i.e., 
A = 256 with A = 2) by locating the maximum of the 
variance of the order parameter fluctuations in the steady 



3 



ioW 




FIG. 1: (color online.) Finite-size scaling of the steady-state 
density auto-correlation function for (Ln,L±) = (54,24), 
(128,32), and (250,40) at the numerically determined "crit- 
ical temperatures" for these three system sizes according to 
Eq. @; the data are averaged over 2000, 550, and 200 re- 
alizations, respectively. The inset depicts the unsealed data. 




FIG. 2: (color online.) Aging scaling plots for the two-time 
density auto-correlation function following a quench to the 
critical point from a high-temperature initial state, on a lat- 
tice with dimensions L\\ = 125000, L ± = 50 (A = A JSLC = 2). 
Left panel: unsealed data vs. t — s; middle and right panels: 
scaling plots, Eq. (|6]), using the JSLC and RDLG exponents, 
respectively. Each curve is averaged over 200 realizations. 



state, resulting in T^ x 



T, 



250x40 



0.773, Tc 28x32 



0.788 3i 



0.782 and 
As shown in Fig.[TJ we thus achieve 



reasonable but certainly not unambiguous data collapse 
using Eq. ([5]) with the JSLC exponents z\\ = | and 
z||C = (d + A - 2 + r?)/(l + A) = § at d = 2 (Table HJ). 
We also performed simulations at the critical tempera- 
ture in the thermodynamic limit T<§° = 1.41 TS q = 0.800 



H I, EH, H El, E3, where T° q = 0.5673 J is the 
equilibrium critical point of the square Ising lattice [i(J , 
and found no noticeable quality difference in the scaling 
collapse. The data in Fig. Q] highlight the difficulty in ob- 
taining accurate critical scaling in the steady state; the 
density correlations display an extended crossover regime 
before reaching their asymptotic behavior [l3| . 

We therefore next venture to determine the critical ex- 
ponents via a careful aging scaling analysis of the den- 
sity auto-correlation relaxation to the non-equilibrium 
steady state. Initiating the simulation runs with a ran- 
dom particle distribution allows us to investigate the out- 
of-equilibrium relaxation regime wherein time translation 
invariance is broken, and two-point correlations become 
explicit functions of two time variables s (referred to as 
waiting time) and t > s [H}. If a system is quenched to 
the critical point, the initial time sheet may in principle 
induce novel infrared singularities, as is the case for the 
relaxational model A with non-conserved order parame- 
ter [l6|. In contrast, for model B and indeed any system 
with a diffusive conserved order parameter field, no new 
divergences emerge, and the critical initial-slip and aging 
regimes are governed by the steady-state scaling expo- 
nents IB - IS I- We may therefore simply add the waiting 



b = s 1 / z then yields the simple aging scaling form of 
the critical density auto-correlation function, 

C(t,8)=8-<f(t/s), ( = (d + A-2 + r 1 )/z . (6) 

Our simulation data for the two-time density auto- 
correlations at are displayed in Fig. [2 for anisotropic 
lattices with A = A JSLC ; simulation runs with A = 
a rdlg yield virtua ii y identical results j^. The left 
panel, with the unsealed data plotted vs. t — s, con- 
firms that time translation invariance is indeed broken; 
in the middle and right-hand panel, respectively, we at- 
tempt to scale the data either using the JSLC or RDLG 
exponent values from Table [l] The JSLC exponents man- 
ifestly yield excellent data collapse, far superior to the 
RDLG values. Moreover, in order to determine the ag- 
ing exponent in an unbiased manner, we assumed simple 
aging and obtained C by means of the minimization tech- 
nique described in Ref. resulting in £ 0.48, in good 
agreement with the JSLC prediction at d = 2. 

Another commonly studied quantity is the anisotropic 
order parameter, defined in the spin representation as 



m(t,T,L±,L\\) 



sin(7r/Lj 
2UT 



L±,L n 

£ 



,k(t) e 



i2irj/Lj 



time s as independent scaling variable to Eq. ([3]); setting 



(7) 

which is sensitive to the density modulations transverse 
to the drive [H, Q . Its finite-size scaling over a range of 
temperatures near the critical point in the steady state is 
known to be sensitive to the precise location of the criti- 
cal temperature and the "critical region" chosen for the 
data analysis [l4[; stationary finite-size scaling attempts 
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FIG. 3: (color online.) Finite-size scaling for the order pa- 
rameter according to Eq. {SJ. The data and Lu x L± symbols 
are the same as in Fig. 0] which depicts the asymptotic tem- 
poral power law. Left panel: excellent data collapse using the 
JSLC exponents; right panel: attempted data collapse with 
RDLG values. Each curve is averaged over 500 realizations. 



have consequently resulted in v arying values o f the or- 
der parameter exponent j3 
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focus on the temporal evolution of the order parameter 
([7]); near T c and beyond mere microscopic time scales, it 
is governed by the universal scaling law 1|| 43| 



(8) 

Identifying b = yields the critical finite-size scaling 



m(t,L,|) = L|| m(t/L^ 



(9) 



We present our Monte Carlo simulation results in Fig. [3J 
employing two-dimensional anisotropic lattices with A = 
A JSLC . We attempted to scale the data via Eq. ([9]) using 
both JSLC and RDLG exponent values, see Table HI the 
JSLC critical exponents Zj| = | and = | clearly 

give superior data collapse in this "initial slip" region. 

To extract the asymptotic temporal power law in the 
initial slip regime, which allows us to investigate very 
large systems, we set b = t x l z w and r = in Eq. ([8]), 



whence m{t,L\\) = m(t x I z L^,t x I z 

stituting Lj_ = ALl /(1+A) with A 



I x ) . Sub- 
const, to becomes a 
function of x = i^^^LV, 1 only. Following the procedure 



outlined in Ref . [17| j we expand the regular scaling func- 
tion rh(x) = mo + fhxx + . . . for id; here mo = as 
alternating occupied/empty initial conditions were cho- 
sen. Ignoring higher-order terms, we arrive at 

m(t, L||)L|| ~t K , k = (1-j8/^|)/z|| . (10) 
Figure 3] demonstrates convincing data collapse with this 
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FIG. 4: (color online.) Scaling plot according to Eq. (|10|) . 
The dashed / solid straight lines above / below the simulation 
data respectively correspond to the predicted RDLG / JSLC 
power laws. The inset shows the unsealed order parameter 
growth. Each curve is averaged over 500 realizations. 



scaling form. Performing a fit to the 64000 x 40 lattice 
data yields k 



0.487, close to the JSLC value § 



44|. 



In summary, we present detailed Monte Carlo simula- 
tion results for the KLS model at criticality confirming 
the JSLC scaling exponents through the measurement of 
various dynamical quantities, specifically the order pa- 
rameter in the initial slip and the two-time density auto- 
correlation function in the aging regime, which can be 
numerically accessed for very large lattices. Finite-size 
scaling of the density auto-correlation function gives sat- 
isfactory (if not fully convincing) data collapse with the 
JSLC exponents. Simple scaling arguments allow us to 
relate the aging exponent to known steady-state criti- 
cal exponents [lg, |43[ . We have performed Monte Carlo 
simulations in the transient regime where time transla- 
tion invariance is broken, and measured the density auto- 
correlation at different waiting times, thus exploring in 
detail the universal relaxation features towards a critical 
non-equilibrium stationary state. The aging exponent in- 
ferred from the JSLC exponent values yields convincing 
data collapse over a range of waiting times, even with dif- 
ferent lattice anisotropies. In contrast, the RDLG scaling 
exponents do not permit a reliable simple aging scaling 
collapse. The time-dependent order parameter data at 
short times for different lattice sizes were found to fit a 
universal master curve using finite-size scaling with JSLC 
exponents. Applying the same scaling form with RDLG 
exponents again could not collapse our data. Supporting 
the conclusions from a careful numerical finite-size scal- 
ing analysis 37J , our work thus provides strong evidence 
that the universal features of the critical KLS model are 
(even with infinite drive) adequately described by the 
JSLC coarse-grained stochastic Langevin equation, and 
captured by the associated field-theoretic RG approach. 
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